\documentclass{article}
\usepackage{geometry}
\geometry{margin=1.5cm, vmargin={0pt,1cm}}
\setlength{\topmargin}{-1cm}
\setlength{\paperheight}{29.7cm}
\setlength{\textheight}{25.3cm}

% useful packages.
\usepackage{amsfonts}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{amsthm}
\usepackage{enumerate}
\usepackage{graphicx}
\usepackage{multicol}
\usepackage{fancyhdr}
\usepackage{layout}
%\usepackage{ctex}
\usepackage{listings}
\usepackage{subfigure}


% some common command
\newcommand{\dif}{\mathrm{d}}
\newcommand{\avg}[1]{\left\langle #1 \right\rangle}
\newcommand{\difFrac}[2]{\frac{\dif #1}{\dif #2}}
\newcommand{\pdfFrac}[2]{\frac{\partial #1}{\partial #2}}
\newcommand{\OFL}{\mathrm{OFL}}
\newcommand{\UFL}{\mathrm{UFL}}
\newcommand{\fl}{\mathrm{fl}}
\newcommand{\op}{\odot}
\newcommand{\Eabs}{E_{\mathrm{abs}}}
\newcommand{\Erel}{E_{\mathrm{rel}}}

\usepackage{xcolor}
\usepackage{fontspec} 
\definecolor{dkgreen}{rgb}{0,0.6,0}
\definecolor{gray}{rgb}{0.5,0.5,0.5}
\definecolor{comment}{rgb}{0.56,0.64,0.68}

\newfontfamily\monaco{Monaco}
\lstset {
aboveskip=3mm,
belowskip=3mm,
showstringspaces=false,       % underline spaces within strings
columns=flexible,
framerule=1pt,
rulecolor=\color{gray!35},
backgroundcolor=\color{gray!5},
basicstyle={\small\monaco},           % the size of the fonts that are used for the code
numbers=left,                   % where to put the line-numbers
numberstyle=\tiny\monaco\color{gray},  % the style that is used for the line-numbers
numbersep=5pt,                  % how far the line-numbers are from the code
commentstyle=\color{comment},
keywordstyle=\color{blue},
stringstyle=\color{dkgreen},
tabsize=2,                      % sets default tabsize to 2 spaces
captionpos=b,                   % sets the caption-position to bottom
breaklines=true,                % sets automatic line breaking
breakatwhitespace=false,        % sets if automatic breaks should only happen at whitespace
escapeinside={\%*}{*)},            % if add LaTeX within your code
morekeywords={*,...}               % if add more keywords to the set
}


\begin{document}
\title{Homework \#4}
\pagestyle{fancy}
\lhead{Name Li HuiTeng 3180102114}
\chead{ numPDE\#4}
\rhead{Date 21.05.20}


\tableofcontents

\newpage

\section{9.8 Prove the relative error is bounded by its relative residual.}
\begin{proof}[Proof]
	For the left side, since
	\[
		\left\|u\right\|\leqslant\left\|A^{-1}\right\|\left\|Au\right\|,\;\left\|Ae\right\|\leqslant\left\|A\right\|\left\|e\right\|,	\]
	we have
	\begin{align*}
		\frac{\vert\vert Ae\vert\vert}{\vert\vert A^{-1}\vert\vert\vert\vert Au\vert\vert} & \leq\frac{\vert\vert A\vert\vert\vert\vert e\vert\vert}{\vert\vert u\vert\vert}, \\\frac1{cond(A)}\frac{\vert\vert r\vert\vert}{\vert\vert f\vert\vert}&=\frac1{cond(A)}\frac{\vert\vert Ae\vert\vert}{\vert\vert Au\vert\vert}\leq\frac{\vert\vert e\vert\vert}{\vert\vert u\vert\vert}.
	\end{align*}
	For the right side, similarly we have
	\begin{align*}
		\frac{\vert\vert Au\vert\vert}{\vert\vert A^{-1}\vert\vert\vert\vert Ae\vert\vert} & \leq\frac{\vert\vert A\vert\vert\vert\vert u\vert\vert}{\vert\vert e\vert\vert}， \\cond(A)\frac{\vert\vert r\vert\vert}{\vert\vert f\vert\vert}&=\vert\vert A\vert\vert\vert\vert A^{-1}\vert\vert\frac{\vert\vert Ae\vert\vert}{\vert\vert Au\vert\vert}\geq\frac{\vert\vert e\vert\vert}{\vert\vert u\vert\vert}.
	\end{align*}
\end{proof}

\section{9.9 Compute the values of cond(A) for n=8,1024.}
\begin{proof}[Proof]
	Since
	\[
		cond(A)=\frac{sin^2\frac{(n-1)\pi}{2n}}{sin^2\frac{\pi}{2n}}=cot^2(\frac{\pi}{2n}),
	\]
	we have
	\[
		n=8:\quad cond(A)=2.5\times10^1,
	\]
	\[
		n=1024:\quad cond(A)=4.2\times10^5.
	\]
\end{proof}

\section{9.12 Plot to show that the range of wavenumbers is [1,n) if the boundary values are 0. If not, show what happens.}
\begin{proof}[Proof]
	WLOG suppose n=10. The figure below has shown an example about how wavenumber k=1 or k=n-1 works. 
	Also, if boundary values are not obliged to be zero, then the range of wavenumbers will be [1,n] with k=n shown below.
	\begin{center}
		\includegraphics[width=0.7\textwidth]{pro912.png}
	\end{center}
	
\end{proof}


\section{9.15 Plot the case of n=6 for Aliasing Lemma.}
\begin{proof}[Proof]
	The figure is shown below.

\begin{center}
	\includegraphics[width=0.7\textwidth]{pro915.png}
\end{center}
\end{proof}


\section{9.19 Prove Iteration Convergence Lemma.}
\begin{proof}[Proof]
	To prove the lemma, we just need to show that $\lim_{k \to \infty} T^k=O$ iff all eigenvalues of
	T satisfy $|\lambda | < 1$.

	By Theorem 7.67, we have
	\[
		T=RJR^{-1}\Rightarrow T^k=RJ^{k}R^{-1}.
	\]
	Therefore $T^k$'s convergence (towards 0) is determined by the convergence (towards 0) of
	\[
		J_i^n=\begin{bmatrix}\lambda_i^n&\begin{pmatrix}n\\1\end{pmatrix}\lambda_i^{n-1}&\begin{pmatrix}n\\2\end{pmatrix}\lambda_i^{n-2}&\cdots&\begin{pmatrix}n\\m_i-1\end{pmatrix}\lambda_i^{n-m_i+1}\\&\lambda_i^n&\begin{pmatrix}n\\1\end{pmatrix}\lambda_i^{n-1}&\cdots&\begin{pmatrix}n\\m_i-2\end{pmatrix}\lambda_i^{n-m_i+2}\\&&\ddots&\ddots&\vdots\\&&&\lambda_i^n&\begin{pmatrix}n\\1\end{pmatrix}\lambda_i^{n-1}\\&&&&\lambda_i^n\end{bmatrix},
	\]
	which follows from $J_i^n=(\lambda_iI+\eta )^n$ where $\eta$ is the nilpotent matrix with $\eta^{m_i}=0$,
	\begin{align*}
		\eta_{ij}=\left\{\begin{array}{lc}1&j-i=1;\\0&otherwise.\end{array}\right. \\
	\end{align*}
	As $k \to \infty$, the nonzero elements of $J_i^n$ approach 0 if and only if $|\lambda_i|<1$, which completes the proof.
\end{proof}

\section{9.24 Verify eigenvalues and eigenvectors of the weighted Jacobi matrix.}
\begin{proof}[Proof]
	For every eigenvector $V_k$ of A, we have
	\[
		T_\omega V_k=(I-{\textstyle\frac{\omega h^2}2}A)V_k=(1-{\textstyle\frac{\omega h^2}2}\lambda_k)V_k.
	\]
	Since $V_k$ are distinct for $k=0,1,\cdots,n-1$ and $dim(T_\omega)=n-1$, we conclude that the eigenvectors of
	$T_\omega$ are the same of those of A, with the corresponding eigenvalues as
	$(1-2\omega sin^2 \frac{k\pi}{2n})$, for $k=0,1,\cdots,n-1$.
\end{proof}

\section{9.25 Reproduce Fig.2.7. and verify the slow convergence for n=64.}
\begin{proof}[Proof]
	The results are reproduced by MATLAB R2019b.
	\lstset{language=Matlab}
	\begin{lstlisting}
%This is a file named pro924.m.
	clear;
	n=64;
	k=0:0.01:n;
	w=[1/3;1/2;2/3;1];
	figure(1);
	for i=1:4
	y=1-2*w(i)*sin(k*pi/(2*n)).^2;
	plot(k,y);
	hold on
	end
	xlabel('k'),ylabel('y'),title('Eigenvalues of the iteration matrix R_ω'),grid on
	saveas(1,'924Eigenvalues.png');
\end{lstlisting}

	\begin{center}
		\includegraphics[width=0.7\textwidth]{eigenvalues924.png}
	\end{center}

	For n=64, $\omega \in [0,1]$, we have
	\[
		\rho(T_\omega)\geq1-2sin^2(\pi/128)=0.9988\geq0.9986
	\]
	Hence, the convergence is quite slow.
\end{proof}

\section{9.28 Reproduce Fig.2.8. and verify the difference between $\omega=1$ and $\omega=2/3$.}
\begin{proof}[Proof]

	The results are produced by MATLAB R2019b.
	\lstset{language=Matlab}
	\begin{lstlisting}
%This is a file named pro927.m.
clear;
clc;
N=64;
A_noh=diag(2*ones(N,1))-diag(1*ones(N-1,1),1)-diag(1*ones(N-1,1),-1);
I=eye(N);
Z=zeros(1,N);
h=1/N;
j=1:N;
k=64;
iteration=zeros(N-1,2);
w=[1;2/3];
a = 1:63;
figure(1);
for l=1:2
Tw=I-(w(l)/2)*A_noh;
for k=1:N-1
omega=sin(j*k*pi*h)';
flag=0;
err_rel=0;
temp1=omega;
omega_norm=norm(omega);
while(err_rel<100 && flag<300)
    temp1=Tw*temp1;
    err_rel=omega_norm/norm(temp1);
   flag=flag+1; 
end
iteration(l,k)=flag;
end
end
b=iteration(1,:);
values = spcrv([[a(1) a a(end)];[b(1) b b(end)]],3);
plot(values(1,:),values(2,:),'r-');
hold on;
b=iteration(2,:);
values = spcrv([[a(1) a a(end)];[b(1) b b(end)]],3);
plot(values(1,:),values(2,:),'b-.');
hold on;
y=0:300;
x=16*ones(1,301);
plot(x,y,'k');
hold on;
x=3*x;
plot(x,y,'c');
hold on;
legend('w=1','w=2/3','x=16','x=48'),xlabel('k'),ylabel('Iterations'),title('Weighted Jacobi method with ω = 1 and ω = 2/3'),grid on;
saveas(1,'927Interations.png');
\end{lstlisting}
	\begin{center}
		\includegraphics[width=0.7\textwidth]{927Interations.png}
		\title{927Interations.png}
	\end{center}
	From the figure above, we conclude that regular Jacobi is only good for damping modes
	$16 \leq k \leq 48$; while for $\omega=2/3$, the modes $16 \leq k \leq 63$ are all damped
	out quickly.
\end{proof}

\section{9.41 Show the computational cost of an FMG cycle is less than $\frac{2}{(1-2^{-D})^2}WU$ and
give the supremum of bounds for D=1,2,3.}
\begin{proof}[Proof]
	As we can see, in a V-cycle the computational cost is 
	\[
		S_V(m)=2WU\sum_{i=0}^{m}2^{-iD}=2WU\frac{1-2^{-(m+1)D}}{1-2^{-D}}.
	\]
	Also, a V-cycle beginning from $\Omega^{2^{i}h}$ costs $2^{−iD}S_V(m-i)$. Since an FMG cycle is composed of single V-cycles 
	beginning from $\Omega^{2^{i}h}$, $i=0,1,\cdots,m-1,m$, the computational cost of an FMG cycle with $n=2^m$ cells can be expressed as
	\[
		S_{FMG}(m)=\sum_{i=0}^{m}2^{-iD}S_V(m-i)\leq S_V(m)\sum_{i=0}^{m}2^{-iD}=2WU\frac{(1-2^{-(m+1)D})^2}{{(1-2^{-D})}^2}\leq\frac{2WU}{{(1-2^{-D})}^2},
	\]
	which completes the proof. The first "$\leq$" is because $S_V$ increases as $m$ increases.

	Next we will determine the upper bounds of $S_{FMG}$ for different dimensions. In fact, we have 
	\begin{align*}
		S_{FMG}(m)&=\sum_{i=0}^{m}2^{-iD}S_V(m-i)\\
		&=\sum_{i=0}^{m}2^{-iD}\cdot 2WU \cdot \frac{1-2^{-(m+1-i)D}}{1-2^{-D}}\\
		&=\frac{2WU}{1-2^{-D}}\sum_{i=0}^{m}(2^{-iD}-2^{-(m+1)D})\\
		&=\frac{2WU}{1-2^{-D}}(\frac{1-2^{-(m+1)D}}{1-2^{-D}}-(m+1)2^{-(m+1)D}).
	\end{align*}

	For D=1:
	\begin{align*}
		S_1=4WU(2-(m+3)2^{-(m+1)}).
	\end{align*}
	Let $a_n=n\cdot2^{-n}$. And we have
	\[
		a_{n+1}<a_n\Leftrightarrow(n+1)2^{-(n+1)}<n2^{-n}\Leftrightarrow1<n.
	\]
	Thus, the infimum of $\{a_n\}_{n \geq 4}$ is $lim_{n\to \infty}a_n=0$. Take $n=m+3$, and the supremum of bound is $8WU$.

	For D=2:
	\begin{align*}
		S_2=\frac{8}{9}WU(4-(3m+7)\cdot4^{-(m+1)}).
	\end{align*}
	Let $a_n=(3n+4)\cdot4^{-n}$. And we have
	\[
		a_{n+1}<a_n\Leftrightarrow(3n+7)4^{-(n+1)}<(3n+4)4^{-n}\Leftrightarrow -1<n.
	\]
	Thus, the infimum of $\{a_n\}_{n \geq 2}$ is $lim_{n\to \infty}a_n=0$ and the supremum of bound is $\frac{16}{9}WU$.

	For D=3:
	\begin{align*}
		S_3=\frac{16}{49}WU(8-(7m+15)\cdot8^{-(m+1)}).
	\end{align*}
	Let $a_n=(7n+8)\cdot8^{-n}$. And we have
	\[
		a_{n+1}<a_n\Leftrightarrow(7n+15)8^{-(n+1)}<(7n+8)8^{-n}\Rightarrow 0<n.
	\]
	Thus, the infimum of $\{a_n\}_{n \geq 1}$ is $lim_{n\to \infty}a_n=0$ and the supremum of bound is $\frac{64}{49}WU$.
\end{proof}

\section{9.47 Explain why the magnitude of all four $c_i's$ are small and reproduce the plots.}
\begin{proof}
	To give a clear and thorough explanation to the question, first let's take $c_1$ as an example to see why $c_1$ is near 0.
	(We only consider the situation that $v_1,v_2$ are non-zero.)

	For $c_1=\lambda_k^{v_1+v_2}s_k$,
	on one side, when k is small, it implies an LF mode. Since Jacobi methods are only good for damping HF modes, $\lambda_k$ should tend
	to 1. However, $s_k=sin^2(\frac{k\pi}{2n})$ is near 0 when k is small and it reduces the norm of $c_1$ greatly.

	On the other side, when
	k is large, it implies an HF mode. Although $s_k$ is near 1, the HF mode, however, will be damped out quite quickly by
	Jacobi methods. Hence, the corresponding eigenvalue $\lambda_k$ has to be rather small and near 0.

	For other three $c_i$, the proof is similar to $c_1's$. In fact, this sort of phenomena should be attributed to a perfect
	match of Jacobi methods and interpolation. Jacobi methods are good at damping HF modes, while interpolation is most effective
	when the error is smooth. Since members of the pair can make up for drawbacks from counterpart, I think they are born to make pair.

	The results are reproduced by MATLAB R2019b.
	\lstset{language=Matlab}
	\begin{lstlisting}
%This is a file named create_plots.m.
	function y=create_plots(v1,v2)
	string=['pro947_',num2str(v1),'_',num2str(v2),'.png'];
	N=64;
	h=1/N;
	w=2/3;
	syms k;
	lam=1-2*w*sin(k*pi/(2*N))*sin(k*pi/(2*N));
	lamc=1-2*w*cos(k*pi/(2*N))*cos(k*pi/(2*N));
	sk=sin(k*pi/(2*N))*sin(k*pi/(2*N));
	ck=cos(k*pi/(2*N))*cos(k*pi/(2*N));
	c1(k)=sk*lam.^(v1+v2);
	c2(k)=sk*(lam.^v1)*(lamc.^v2);
	c3(k)=ck*(lam.^v2)*(lamc.^v1);
	c4(k)=ck*lamc.^(v1+v2);
	figure(1);
	fp1 = fplot(@(k)c1(k),[0 32]);
	grid on;
	fp1.LineStyle = '-';
	fp1.Color = 'k';
	hold on;
	fp2 = fplot(@(k)c2(64-k),[32 64]);
	fp2.Color = 'b';
	fp2.LineStyle = '-.';
	fp3 = fplot(@(k)c3(k),[0,32]);
	fp3.Color = 'm';
	fp3.LineStyle = '--';
	fp4 = fplot(@(k)c4(64-k),[32 64]);
	fp4.Color = 'r';
	fp4.LineStyle = '-';
	fp4.Marker = '.';
	title(string);
	legend('c1','c2','c3','c4');	
	saveas(1,string);
	clf(1);
	y=0;
	end

%This is a file named pro947.m.
a(1,:)=[0,0,1,2,2,4];
a(2,:)=[0,2,1,0,2,0];
for i=1:6
create_plots(a(1,i),a(2,i));
end
\end{lstlisting}
	The output is shown below.
	\newpage
	\includegraphics[width=0.45\textwidth]{pro947_0_0.png}
	\includegraphics[width=0.45\textwidth]{pro947_0_2.png}

	\includegraphics[width=0.45\textwidth]{pro947_1_1.png}
	\includegraphics[width=0.45\textwidth]{pro947_2_0.png}

	\includegraphics[width=0.45\textwidth]{pro947_2_2.png}
	\includegraphics[width=0.45\textwidth]{pro947_4_0.png}
\end{proof}

\section{9.51 Prove $\dim \mathcal{R}(I^{2h}_{h})=\frac{n-2}{2}$, $\dim \mathcal{N}(I^{2h}_{h})=\frac{n}{2}$. }
\begin{proof}[Proof]
	This follows directly from (9.39) and Lemma 9.49. Applying fundamental theorem of linear algebra to above results, we have
	\begin{align*}
		\dim \mathcal{R}(I^{2h}_{h}) & =\frac{n-2}{2}-\dim \mathcal{N}({(I^{2h}_{h})}^T)=\frac{n-2}{2}-\dim \mathcal{N}(I^{h}_{2h})=\frac{n-2}{2}, \\
		\dim \mathcal{N}(I^{2h}_{h}) & =\dim \mathcal{N}({(I^{h}_{2h})}^T)=n-1-\dim \mathcal{R}(I^{h}_{2h})=n-1-\frac{n-2}{2}=\frac{n}{2}.
	\end{align*}
\end{proof}

\section{9.56 Show (9.45) holds.}
\begin{proof}[Proof]
	This follows directly from (9.39). We have
	\begin{align*}
		||v^{2h}-u^{2h}||_{A^{2h}} & =(A^{2h}(v^{2h}-u^{2h}),v^{2h}-u^{2h})                       \\
		                           & =(I^{2h}_{h}A^{h}I^{h}_{2h}(v^{2h}-u^{2h}),v^{2h}-u^{2h})    \\
		                           & =c(A^{h}I^{h}_{2h}(v^{2h}-u^{2h}),I^{h}_{2h}(v^{2h}-u^{2h})) \\
		                           & =c||I^{h}_{2h}(v^{h}-u^{h})||_{A^h}.
	\end{align*}
\end{proof}

\end{document}

%%% Local Variables: 
%%% mode: latex
%%% TeX-master: t
%%% End: 
